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The Navier-Stokes equations and boundary conditions for viscous fluids of capillary size are for- 
mulated in curvilinear coordinates associated with a geometry of the fiuid-gas interface. As a result, 
the fluid dynamics of drops and menisci can be described taking into account an influence of gravi- 
tational forces and flows on the surface shape. This gives a convenient basis for respective numerical 
studies. Estimations of the effects are presented for the case of an evaporating sessile drop. 



Introduction. A number of important physical features in studying fluid flows in evaporating liquid drops and 
menisci of capillary size has been found recently both theoretically and experimentally [Ij, J,, 0, H SB- In 
particular, it was demonstrated that the vortex convection takes place in evaporating drops and menisci under various 
conditions The activity in the field is associated now with important applications. The particular examples 

are the evaporative contact line deposition P,[3,f4, 9, 10, 11], studies of DNA stretching behavior and DNA mapping 
methods fl3. ITsl. 13 . developing methods for jet ink printing [T^, [l^ self-assembly of nanocrystal superlattice 
monolayerlialMHoi • 

For describing the processes theoretically one should carry out, in general, a joint study of the fluid dynamics, 
the thermal conduction and the vapor diffusion together with respective boundary conditions, in particular at the 
fluid-gas interface. Standard approximations used in the theoretical studies are a spherical cap shape of the drop or 
menisci and a neglection of the hydrodynamical pressure and velocities in the generalized Laplace formula. Though 
such approximations can be justified under certain conditions, there are a wide range of parameters of the problem 
when a more accurate theoretical description of liquid surfaces of capillary size is needed. 

A shape of a surface is, generally, controlled by combined effects of a surface tension, gravitational forces, a 
hydrodynamic pressure and a velocity distribution near the surface. For solving fluid dynamics problems in a vicinity 
of curved surfaces of an arbitrary shape, an explicit approach is developed in the present paper, making use of "natural" 
curvilinear coordinates associated with a surface geometry. Both fluid dynamics equations and the respective boundary 
conditions are formulated in these coordinates. The equations in such a form are convenient for numerical simulations. 
We also present analytical estimations for the effects in question for the case of an evaporating sessile drop, which 
follow from the obtained results. 

Equations and boundary conditions. The Navier-Stokes equations take the form 

dv 1 

— — h (vV)v H — gradp = i/Av. (1) 
ot p 

For simplicity, we assume below that a shape of the surface as well as fluid flows are axially symmetric and vg — 0, 
where r, 9, z are cylindrical coordinates. This property is valid for a wide class of problems. Therefore, it is 
convenient to use cylindrical coordinates in the bulk of an incompressible viscous liquid and to introduce vorticity 
7 = dvr/dz — dvz/dr and the stream function -0, such that d^/dz = rvr, d^/dr — —rvz (as distinct from the 
two-dimensional case [2l|b Then rot v ~ 7(r, z)ig, the continuity equation divv = is naturally satisfied. Equations 
for quantities 7, tp are given by 

^7(r, z) + (vV)7(r, z) = (^A7(r, z) - ^^^^ , (2) 

A^) — = r7. (3) 

r Or 

In order to formulate equations close to the surface, it is convenient to choose orthogonal curvilinear coordinates 
x"{x,y,z), x'^{x,y,z), x^{x,y,z) with local basis vectors normal and tangential to the surface at every point. In 
order to write down in these curvilinear coordinates the differential forms which enter the hydrodynamic equations, 
one needs to find explicit expressions for the metric tensor and Christoffel symbols for the chosen class of coordinate 
systems. Consider both the contravariant coordinates and the respective physical curvilinear coordinates 

n, r, Tg. Locally dn is a length along the normal to the surface, dr is a surface arc length in the meridian plane, and 
drg is a surface arc length associated with the rotation angle around the z axis. For an axially symmetric surface 
drg = r(n, T)dO. For a differential of radius-vector we have 

dr = dx^Gn + dx'^e-r + dx^eg = dnin + dr -I- rdd ig , (4) 
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where = idx/dx^ +idy/dx^ + kdz/dx^ are contravariant base vectors. Unlike contravariant base vectors e„,eT-,ee, 
Cartesian base vectors i,j,k and physical curvilinear base vectors i„, i,-, ig are orthonormalized i„ = e„/ y/gnn, 
ir = er/ ^/grr, 10 = Gg/y/ggg. The validity of the following relations is necessary for constucting the contravariant 
basis 



dej 
dx^ 



dej_ 

dx^ 



(5) 



The unit vectors of physical coordinate system do not satisfy such a requirement, in contrast to contravariant basis 
vectors, due to the difference in their normalizations. We note that the requirement ([5]) will be satisfied if the local 
angle ip between the normal vector to the surface and the symmetry axis depends only on x'^ , doesn't depend on x" 
and x^, and 



dr 

dx'^ 



dip 



- + /o(x^), i.e 



dip 



(6) 



dx'^ ' -"j^ ~ dx'^dn dx"^ 

Here function fo{x'^) is defined by the geometry of the problem and, in particular, by a choice of the origin for the 
coordinate n. For a spherical surface x"^ = (p, x" = n = R = \J x^ -V -V ^ x'^ = ip = t / n, x^ = and r = nsunp. 
Relations ([6]) permit to determine the contravariant basis near the surface: e„ = i„, = irdr / dx'^ , eg = r\g. 



One obtains for the basis the following components of the metric tensor g^ 



1, grT = (dr/dx^) , ggg 



0, g = det gik = r [dr/dx'^) , and the corresponding Christoffel symbols 



dx-^dr/dx^' 



ru T^" 



smtp 



dx^ dx^ ' 



d'r 1 



-L — i At — 



cos(^ 9t 
r dx'^ ' 



Tee = sin ^3, 



r flfl — 



r cos tp 
dr I dx'^ 



(7) 



The expressions for the metric tensor and Christoffel symbols allow to obtain explicit formulas for all differential 
forms, according to general rules of the differential geometry [2^ . In particular, for arbitrary vector F one finds 



rotr 



1 fd{rFg) dFr 



86 



1 /9F„ d{rFg) 



\ do 



dr 



dn 



dtp 



dr dr 



Fr ig 



(8) 



Therefore, 



7 
Av 

(vV)v 



(rotv)e = 

-rot (jie) 
dvn 

'Jn-^ h 

on 



dVr _ dVn 

dn dr 



dip 



57 + ^^'+^^ 



97 cos tp 

dvn dip 
dr dr 



sm p 



dj 
dn 
dvr 
dn 



Thus, the components of Eq.([T]) may be rewritten as 



dp 1 


^dVr 


dVr 


d^ = -P\ 


V dt 


dr 


dp 


( dvn 






V dt 


dn 



dVn 

dr 

dVr 

dn 



+ 7 

- 7 



97 
d-f 



dVr 
57 



sm(^ 
r 

cos ip 



dip 
'd^ 



(9) 
(10) 

(11) 

(12) 
(13) 



In a more general case when vg ^ 0, the terms pvg cos ip/r, pvgsinip/r should be added to right-hand member of 
Eqs. p^ . P^ correspondingly. 

The components of a viscous stress tensor ct^^ — r]{dvi / dxk + dvk / dxt) , which describe momentum transfer through 
the boundary, take the form 



dn 



dVn ^ dVr dip 

dr dn ^ dr 



(14) 



The boundary condition at the surface is |2l| 



P-Pv-<J[^^ 



1 \\ f , 9vk\ da 

I I ni = T] I \ 1 



i?2 



V dxk dxi 



dxi ' 



(15) 
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The normal vector is directed here in the outward direction, towards the atmosphere. Here py a is pressure of the gas 
and atmosphere, P is a hydrodynamic pressure on the surface, Ri^2 are main local radii of curvature on the surface, 
and (T is a surface tension. Projections of (fT5|) to the local tangential and normal directions to the surface are 



da 

P-Pv = 



dVn dVr dip 

dr dn ^ rfr 



1 1 A „ 9u„ 



(16) 
(17) 



Eg. pop is the boundary condition for the velocities on the surface. In particular, one gets the boundary condition for 
the vorticity at the surface: 



1 da „ 
7^ -— + 2 
rj dr 



dip dvn 
^ dr dr 



(18) 



Boundary conditions for stream function at the surface may be obtained by integrating the expression dip/dr — 
—rVn{T), where Vnir) is the normal component of the velocity to the boundary. The boundary conditions for the 
stream function are particulary simple if the motion of the surface is much slower than typical fluid velocities of a 
problem, when one can put w„ w 0. 

Eq. pT|) represents the boundary condition that allows to obtain a shape of the surface. The pressure P{t) satisfies 
Navier-Stokes equations p^ . pS]) with corresponding projections of the gravitational force added to the right-hand 
part of the equations. Therefore, the quantity p{t) — P{t) + pgz satisfies Egs. p^ . (fO)l without the additional 
terms. To further simplify equations and boundary conditions, one can introduce the quantities P3(t) — P3(0) 



p(r) 

ko = 



p{Q) - 27]dvn/dn\l, p4(r) 



■A. 



P4(0) = p{t) 
Then we get 



p(0) + p{vl -|-u2)/2|p, and k = R-^^ + R^^ = dip/dr + s\nip/r, 



Pair) 
dp 

d^ ~ 



ko 



P3(0) ^a{k~ko)+ pg{z - zq), 
P3(t) -P3(0) - pg{z - Zo) 



smp 



(19) 
(20) 



In the particular case when one can disregard the term with the pressure, Eq. ([20|) reduces to the Young-Laplace 
equation in the form obtained in ^23i] . The tangential component (|12p of the Navier-Stokes equation may be represented 
as 



dpi 
dr 



-P 



dVr 

~dt 



V 



dn 



sm p) 
-TT- H 7 



(21) 



The equation (PO)) turns out to be quite convenient for determining the shape of the surface. The shape of an axially 
symmetric surface is unambiguously described by the function (/^(t). Because all the expressions contain either the 
diff'erence P4(t) — p4(0) or the derivative dpi{T)/dT, an initial value of Pi{Q) is still an arbitrary constant. It is 
convenient to take 



pM = 277 



dn 



T=0 



T = 



Then 



Pz{t) ~ P3(0) = Pa{t) - ^ ■ 



(22) 



(23) 



Introducing the vector y = (r(T), (/^(t), z(t),P4(t))^ allows to represent Eqs. (^0)1 . (PT|) . . dr(r)/dT — cosip, 
dz{T)/dT ~ — smp) in the following form 



ft \ 

= f{r,y)- 



(24) 



At the initial point one has y(0) = (r(0), (^(0), z(0),p4(0)) . Here p4(0) is defined in (HH). The Cauchy problem for 
the system of differential equations (|24p with initial conditions derived above can be solved by standard numerical 
methods to obtain the surface profile. 

Estimations. It is of interest to find out a relative role of terms in Eq. pop under specific physical conditions. Below 
we carry out the respective estimatons for an evaporating sessile drop lying on a substrate in the regime of a pinned 
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contact line. Evaporation results in an inhomogeneous spatial temperature distribution in the drop and along the 
drop surface. The corresponding Marangoni forces result in vortex flows of the liquid in the drop. 

The motion of the surface is considered to be much slower than typical fluid velocities. This property is valid 
for a wide class of evaporating drops. Then one can take approximately ?;„ « 0. The fluid motion is considered 
as a quasistationary vortex flow. In the following expressions hq is the characteristic distance between the surface 
of the drop and the vortex center, ro is the contact line radius, a' = —da/dT, AT is the temperature difference 
between the substrate and the apex of the drop, 9c is the contact angle. Therefore, dLp/dr « sai9c/rQ, da/{r]dT) w 
—a'lS.TsvsY9c/{'qrQ9c)- Using the condition (jl6p and taking riQ as a characteristic distance for a change of Vr along 
the normal to the surface, one obtains Vr ~ {dvT-/dn)no — noda / (jidr) + noVrdip/dr, i.e. Vr (1 — nosin^c/fo) ~ 
— a'noAT sin 9c/ {rirQ9c), hence 



< 



a'ATno 
Wo 



(25) 



Therefore \vrdip/dT\ w \vr\s\n9c/ro w nocr'ATsin^ O c / [wl^^ c) < a' AT / (rjro) , i.e. the term Vrdip/dr in (HH) and ([18 
is much smaller than da/{r]dT). It follows from p8|) that 



It follows from \d'^Vr/dn^\ w |fr|/'T-o ^^'^ 



a'AT 
ivro) ' 



d^Vr dip da 



dn dn"^ dr -qdr 



and 1251) that 



dn 



a'AT /sin 



no sin 9^ 
2ro 



a'AT 
Wono ' 



(26) 



(27) 



(28) 



We substitute (pS)) and to (PT|) and integrate the obtained expression over r. This gives the estimation of relative 
effects of pressures and velocities as compared with gravitational forces in Eq. 



Pi{T) -P4(0) 

\pi{r)-pm\ 
pgh 



a'AT9c 
no sin 9c ' 

a'AT 9c 
pgnoh sin 6'c ' 



27? — 

dn 



<P4(t) -P4(0), 



\pvl/2\ ^ 1 [no a' AT 



pgh ^ 2gh \ 77 ro 



(29) 
(30) 



For estimating the term pv^/2 we used (^5]). 

The ratio of gravitational force to the term with surface tension in p9)) is characterised by dimensionless number 
Bo — pghrol {2a suy9c), which is analogous to Bond number. Therefore, pop may be represented as 



\pi{T) - pi{0)\ <j'AT9c ro 



\a(k~ko)\ 



2a sin^ 9c no ' 



\a{k ~ ko)\ 4ro(Tsin0c 



noa'AT 
V 



(31) 



Conclusion. Based on a geometry of the fluid surface, we have derived Eqs. (jl2p .(|13 p and the boundary condition 
(|19p . which allow to obtain numerically a surface profile dynamics and to take into account the influence of fluid 
dynamics and gravitational forces on the shape of the fluid-gas interface. 

The equations and boundary conditions derived in this paper were used in [24| to find numerically the profile of the 
evaporating sessile drop surface. According to Eq. (|30p . the effects of the pressure become more important with the 
increase of the temperature drop in the liquid and with the temperature derivative of the surface tension. Analytical 
estimations ([50)1 applied to the conditions of [23| show that the relative contribution of pressures and velocities as 
compared with gravitational forces in the Laplace formula ([T^ . is not too large. Numerical results confirm this 
qualitative conclusion and give approximately one tenth for the value of this quantity under the conditions of [2J|. 

The author is grateful to V.V. Lebedev and L.N. Shchur for useful discussions and remarks. 
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